subroutine alkane_init( moltype )
use prec
use cell
use general
implicit none
   !----------------------------------------------------------------------------
   integer, intent(in) :: moltype
   !----------------------------------------------------------------------------
   integer  :: i
   !----------------------------------------------------------------------------
   subnamep = subname
   subname  = 'alkane_init'
   !
   if ( allocated( bondsMol  ) ) deallocate( bondsMol  )
   if ( allocated( anglsMol  ) ) deallocate( anglsMol  )
   if ( allocated( dihrsMol  ) ) deallocate( dihrsMol  )
   if ( allocated( attypMol  ) ) deallocate( attypMol  )
   if ( allocated( atposMol  ) ) deallocate( atposMol  )
   if ( allocated( molIDMol  ) ) deallocate( molIDMol  )
   if ( allocated( rotatedMol) ) deallocate( rotatedMol)
   if ( allocated( chargeMol ) ) deallocate( chargeMol )
   if ( allocated( EMass     ) ) deallocate( EMass     )
   if ( allocated( sigma     ) ) deallocate( sigma     )
   if ( allocated( epison    ) ) deallocate( epison    )
   if ( allocated( bondcoef  ) ) deallocate( bondcoef  )
   if ( allocated( anglcoef  ) ) deallocate( anglcoef  )
   if ( allocated( torscoef  ) ) deallocate( torscoef  )
   if ( allocated( DistMin   ) ) deallocate( DistMin   )
   if ( allocated( EName     ) ) deallocate( EName     )
   !
   ntypeMol   = 2  ! 1: CT, Ct; 2: HT, Ht
   nbndtypMol = 2  ! 1: CT-CT;  2: CT-HT
   nangtypMol = 3  ! 1: CT-CT-CT; 2: CT-CT-HT; 3: HT-CT-HT;
   ndihtypMol = 2  ! 1: CT-CT-CT-CT; 2: CT-CT-CT-HT, HT-CT-CT-HT
   !
   allocate( EMass(ntypeMol), sigma(ntypeMol), epison(ntypeMol),           &
   &         bondcoef(2, nbndtypMol), anglcoef(2, nangtypMol),             &
   &         torscoef(4, ndihtypMol), DistMin(ntypeMol), EName(ntypeMol+1) )
   EMass( 1 ) = 12.0107D0
   EMass( 2 ) = 1.00794D0
   sigma      = (/ 3.50D0,  2.50D0 /)
   epison     = (/ 0.066D0, 0.030D0 /)
   DistMin    = (/ 1.7D0, 1.1D0 /)
   EName      = (/ 'C ', 'H ',  'Al' /)
   !
   bondcoef = reshape( (/ 268.D0, 1.529D0,  & ! 1 CT-CT
   &                      340.D0, 1.09D0/), & ! 2 CT-HT
   &                     (/2, nbndtypMol/) )
   anglcoef = reshape( (/ 58.35D0,112.7D0,  & ! 1 CT-CT-CT
   &                      37.5D0, 110.7D0,  & ! 2 CT-CT-HT
   &                      33.D0,  107.8D0/),& ! 3 HT-CT-HT
   &                      (/2, nangtypMol/) )
   torscoef = reshape( (/ 1.300D0, -0.050D0,  0.200D0,  0.0D0, & ! 1 CT-CT-CT-CT CT-CT-CT-CM
   &                      0.000D0,  0.000D0,  0.300D0,  0.0D0/)& ! 2 HT-CT-CT-CT HC-CT-CT-HC
   &                      , (/4, ndihtypMol/) )
   !
   select case ( moltype )
   case ( 6 ) ! Hexane: H3C-CH2-CH2-CH2-CH2-CH3
      !                   1 2   3   4   5   6
      !=========================================================================
      molform    = 'H3C-CH2-CH2-CH2-CH2-CH3'
      molname    = 'Hexane'
      den0       = 0.6548D0
      natomMol   = 20
      nbondsMol  = 19
      nanglsMol  = 36
      ndihrsMol  = 45
      !
      allocate( bondsMol(3, nbondsMol), AnglsMol(4, nanglsMol),                &
      &         dihrsMol(5, ndihrsMol), attypMol(natomMol),                    &
      &         atposMol(3,natomMol), chargeMol(natomMol), molIDMol(natomMol) )                
      !
      attypMol = (/ 1, 2, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 2 /)
      atposMol(:, 1) = (/   -9.9600D0,  1.8300D0,  0.0000D0 /)
      atposMol(:, 2) = (/  -10.7987D0,  2.5261D0,  0.0000D0 /)
      atposMol(:, 3) = (/  -10.0094D0,  1.2024D0,  0.8899D0 /)
      atposMol(:, 4) = (/  -10.0094D0,  1.2024D0, -0.8899D0 /)
      atposMol(:, 5) = (/   -8.6373D0,  2.6187D0,  0.0000D0 /)
      atposMol(:, 6) = (/   -8.6462D0,  3.2477D0, -0.8901D0 /)
      atposMol(:, 7) = (/   -8.6462D0,  3.2477D0,  0.8901D0 /)
      atposMol(:, 8) = (/   -7.2930D0,  1.8676D0,  0.0000D0 /)
      atposMol(:, 9) = (/   -7.2841D0,  1.2385D0,  0.8901D0 /)
      atposMol(:,10) = (/   -7.2841D0,  1.2385D0, -0.8901D0 /)
      atposMol(:,11) = (/   -5.9703D0,  2.6563D0,  0.0000D0 /)
      atposMol(:,12) = (/   -5.9792D0,  3.2853D0, -0.8901D0 /)
      atposMol(:,13) = (/   -5.9792D0,  3.2853D0,  0.8901D0 /)
      atposMol(:,14) = (/   -4.6260D0,  1.9051D0,  0.0000D0 /)
      atposMol(:,15) = (/   -4.6171D0,  1.2761D0,  0.8901D0 /)
      atposMol(:,16) = (/   -4.6171D0,  1.2761D0, -0.8901D0 /)
      atposMol(:,17) = (/   -3.3033D0,  2.6938D0,  0.0000D0 /)
      atposMol(:,18) = (/   -2.4646D0,  1.9977D0,  0.0000D0 /)
      atposMol(:,19) = (/   -3.2539D0,  3.3214D0,  0.8899D0 /)
      atposMol(:,20) = (/   -3.2539D0,  3.3214D0, -0.8899D0 /)
      !
      chargeMol = (/-0.180D0,  0.060D0,  0.060D0,  0.060D0, -0.120D0,  &
      &              0.060D0,  0.060D0, -0.120D0,  0.060D0,  0.060D0,  &
      &             -0.120D0,  0.060D0,  0.060D0, -0.120D0,  0.060D0,  &
      &              0.060D0, -0.180D0,  0.060D0,  0.060D0,  0.060D0  /) 
      !                  Type,  atom 1, atom 2
      bondsMol(:, 1) = (/  2, 2, 1 /)
      bondsMol(:, 2) = (/  2, 3, 1 /)
      bondsMol(:, 3) = (/  2, 4, 1 /)
      bondsMol(:, 4) = (/  1, 5, 1 /)
      bondsMol(:, 5) = (/  2, 6, 5 /)
      bondsMol(:, 6) = (/  2, 7, 5 /)
      bondsMol(:, 7) = (/  1, 8, 5 /)
      bondsMol(:, 8) = (/  2, 9, 8 /)
      bondsMol(:, 9) = (/  2,10, 8 /)
      bondsMol(:,10) = (/  1,11, 8 /)
      bondsMol(:,11) = (/  2,12,11 /)
      bondsMol(:,12) = (/  2,13,11 /)
      bondsMol(:,13) = (/  1,14,11 /)
      bondsMol(:,14) = (/  2,15,14 /)
      bondsMol(:,15) = (/  2,16,14 /)
      bondsMol(:,16) = (/  1,17,14 /)
      bondsMol(:,17) = (/  2,18,17 /)
      bondsMol(:,18) = (/  2,19,17 /)
      bondsMol(:,19) = (/  2,20,17 /)
      !                   Type,  atom1,  atom2, atom3
      AnglsMol(:, 1) = (/ 3, 2, 1, 3 /)
      AnglsMol(:, 2) = (/ 3, 2, 1, 4 /)
      AnglsMol(:, 3) = (/ 2, 2, 1, 5 /)
      AnglsMol(:, 4) = (/ 3, 3, 1, 4 /)
      AnglsMol(:, 5) = (/ 2, 3, 1, 5 /)
      AnglsMol(:, 6) = (/ 2, 4, 1, 5 /)
      AnglsMol(:, 7) = (/ 2, 1, 5, 6 /)
      AnglsMol(:, 8) = (/ 2, 1, 5, 7 /)
      AnglsMol(:, 9) = (/ 1, 1, 5, 8 /)
      AnglsMol(:,10) = (/ 3, 6, 5, 7 /)
      AnglsMol(:,11) = (/ 2, 6, 5, 8 /)
      AnglsMol(:,12) = (/ 2, 7, 5, 8 /)
      AnglsMol(:,13) = (/ 2, 5, 8, 9 /)
      AnglsMol(:,14) = (/ 2, 5, 8,10 /)
      AnglsMol(:,15) = (/ 1, 5, 8,11 /)
      AnglsMol(:,16) = (/ 3, 9, 8,10 /)
      AnglsMol(:,17) = (/ 2, 9, 8,11 /)
      AnglsMol(:,18) = (/ 2,10, 8,11 /)
      AnglsMol(:,19) = (/ 2, 8,11,12 /)
      AnglsMol(:,20) = (/ 2, 8,11,13 /)
      AnglsMol(:,21) = (/ 1, 8,11,14 /)
      AnglsMol(:,22) = (/ 3,12,11,13 /)
      AnglsMol(:,23) = (/ 2,12,11,14 /)
      AnglsMol(:,24) = (/ 2,13,11,14 /)
      AnglsMol(:,25) = (/ 2,11,14,15 /)
      AnglsMol(:,26) = (/ 2,11,14,16 /)
      AnglsMol(:,27) = (/ 1,11,14,17 /)
      AnglsMol(:,28) = (/ 3,15,14,16 /)
      AnglsMol(:,29) = (/ 2,15,14,17 /)
      AnglsMol(:,30) = (/ 2,16,14,17 /)
      AnglsMol(:,31) = (/ 2,14,17,18 /)
      AnglsMol(:,32) = (/ 2,14,17,19 /)
      AnglsMol(:,33) = (/ 2,14,17,20 /)
      AnglsMol(:,34) = (/ 3,18,17,19 /)
      AnglsMol(:,35) = (/ 3,18,17,20 /)
      AnglsMol(:,36) = (/ 3,19,17,20 /)
      !                      Type, atom1, atom2, atom3, atom4
      DihrsMol(:, 1) = (/   2, 6, 5, 1, 2 /)
      DihrsMol(:, 2) = (/   2, 6, 5, 1, 3 /)
      DihrsMol(:, 3) = (/   2, 6, 5, 1, 4 /)
      DihrsMol(:, 4) = (/   2, 7, 5, 1, 2 /)
      DihrsMol(:, 5) = (/   2, 7, 5, 1, 3 /)
      DihrsMol(:, 6) = (/   2, 7, 5, 1, 4 /)
      DihrsMol(:, 7) = (/   2, 8, 5, 1, 2 /)
      DihrsMol(:, 8) = (/   2, 8, 5, 1, 3 /)
      DihrsMol(:, 9) = (/   2, 8, 5, 1, 4 /)
      DihrsMol(:,10) = (/   2, 9, 8, 5, 1 /)
      DihrsMol(:,11) = (/   2, 9, 8, 5, 6 /)
      DihrsMol(:,12) = (/   2, 9, 8, 5, 7 /)
      DihrsMol(:,13) = (/   2,10, 8, 5, 1 /)
      DihrsMol(:,14) = (/   2,10, 8, 5, 6 /)
      DihrsMol(:,15) = (/   2,10, 8, 5, 7 /)
      DihrsMol(:,16) = (/   1,11, 8, 5, 1 /)
      DihrsMol(:,17) = (/   2,11, 8, 5, 6 /)
      DihrsMol(:,18) = (/   2,11, 8, 5, 7 /)
      DihrsMol(:,19) = (/   2,12,11, 8, 5 /)
      DihrsMol(:,20) = (/   2,12,11, 8, 9 /)
      DihrsMol(:,21) = (/   2,12,11, 8,10 /)
      DihrsMol(:,22) = (/   2,13,11, 8, 5 /)
      DihrsMol(:,23) = (/   2,13,11, 8, 9 /)
      DihrsMol(:,24) = (/   2,13,11, 8,10 /)
      DihrsMol(:,25) = (/   1,14,11, 8, 5 /)
      DihrsMol(:,26) = (/   2,14,11, 8, 9 /)
      DihrsMol(:,27) = (/   2,14,11, 8,10 /)
      DihrsMol(:,28) = (/   2,15,14,11, 8 /)
      DihrsMol(:,29) = (/   2,15,14,11,12 /)
      DihrsMol(:,30) = (/   2,15,14,11,13 /)
      DihrsMol(:,31) = (/   2,16,14,11, 8 /)
      DihrsMol(:,32) = (/   2,16,14,11,12 /)
      DihrsMol(:,33) = (/   2,16,14,11,13 /)
      DihrsMol(:,34) = (/   1,17,14,11, 8 /)
      DihrsMol(:,35) = (/   2,17,14,11,12 /)
      DihrsMol(:,36) = (/   2,17,14,11,13 /)
      DihrsMol(:,37) = (/   2,18,17,14,11 /)
      DihrsMol(:,38) = (/   2,18,17,14,15 /)
      DihrsMol(:,39) = (/   2,18,17,14,16 /)
      DihrsMol(:,40) = (/   2,19,17,14,11 /)
      DihrsMol(:,41) = (/   2,19,17,14,15 /)
      DihrsMol(:,42) = (/   2,19,17,14,16 /)
      DihrsMol(:,43) = (/   2,20,17,14,11 /)
      DihrsMol(:,44) = (/   2,20,17,14,15 /)
      DihrsMol(:,45) = (/   2,20,17,14,16 /)
      !
   case ( 12 ) ! hydrogenated hexene dimer
      !                 CH3
      !                 |
      ! CH3-CH2-CH2-CH2-CH-CH-CH2-CH2-CH2-CH3
      !                    |
      !                    CH3
      !=========================================================================
      molform    = 'CH3-(CH2)3-CHCH3-CHCH3-(CH2)3-CH3'
      molname    = 'HexeneDimer'
      den0       = 0.8D0
      natomMol   = 38
      nbondsMol  = 37
      nanglsMol  = 72
      ndihrsMol  = 99
      !
      allocate( bondsMol(3, nbondsMol), AnglsMol(4, nanglsMol),                &
      &         dihrsMol(5, ndihrsMol), attypMol(natomMol),                    &
      &         atposMol(3,natomMol), chargeMol(natomMol), molIDMol(natomMol) )
      !
      attypMol = (/ 1, 2, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 1, 2, 1, 2, 2,&
      &             1, 2, 2, 1, 2, 2, 1, 2, 2, 2, 1, 2, 2, 2, 1, 2, 2, 2 /)
      atposMol(:, 1) = (/  -9.5596D0,  1.2736D0, 0.0799D0 /)
      atposMol(:, 2) = (/ -10.4145D0,  1.8734D0,-0.2654D0 /)
      atposMol(:, 3) = (/  -9.6317D0,  1.1231D0, 1.1667D0 /)
      atposMol(:, 4) = (/  -9.5671D0,  0.2978D0,-0.4269D0 /)
      atposMol(:, 5) = (/  -8.2520D0,  2.0095D0,-0.2480D0 /)
      atposMol(:, 6) = (/  -8.1762D0,  2.1601D0,-1.3345D0 /)
      atposMol(:, 7) = (/  -8.2418D0,  2.9850D0, 0.2596D0 /)
      atposMol(:, 8) = (/  -7.0644D0,  1.1696D0, 0.2357D0 /)
      atposMol(:, 9) = (/  -7.1395D0,  1.0220D0, 1.3223D0 /)
      atposMol(:,10) = (/  -7.0878D0,  0.2073D0,-0.2689D0 /)
      atposMol(:,11) = (/  -5.7568D0,  1.9270D0,-0.0979D0 /)
      atposMol(:,12) = (/  -5.7202D0,  2.1001D0,-1.1827D0 /)
      atposMol(:,13) = (/  -5.8367D0,  2.8966D0, 0.4153D0 /)
      atposMol(:,14) = (/  -4.3552D0,  1.2860D0, 0.3525D0 /)
      atposMol(:,15) = (/  -4.3542D0,  1.1616D0, 1.4448D0 /)
      atposMol(:,16) = (/  -3.9166D0, -0.1767D0,-0.3525D0 /)
      atposMol(:,17) = (/  -3.9176D0, -0.0522D0,-1.4447D0 /)
      atposMol(:,18) = (/  -2.5150D0, -0.8177D0, 0.0980D0 /)
      atposMol(:,19) = (/  -2.5516D0, -0.9908D0, 1.1827D0 /)
      atposMol(:,20) = (/  -2.4351D0, -1.7872D0,-0.4153D0 /)
      atposMol(:,21) = (/  -1.2074D0, -0.0602D0,-0.2357D0 /)
      atposMol(:,22) = (/  -1.1323D0,  0.0873D0,-1.3222D0 /)
      atposMol(:,23) = (/  -1.1840D0,  0.9021D0, 0.2689D0 /)
      atposMol(:,24) = (/  -0.0198D0, -0.9001D0, 0.2480D0 /)
      atposMol(:,25) = (/  -0.0956D0, -1.0507D0, 1.3346D0 /)
      atposMol(:,26) = (/  -0.0300D0, -1.8756D0,-0.2595D0 /)
      atposMol(:,27) = (/   1.2878D0, -0.1642D0,-0.0799D0 /)
      atposMol(:,28) = (/   2.1427D0, -0.7640D0, 0.2654D0 /)
      atposMol(:,29) = (/   1.2953D0,  0.8115D0, 0.4269D0 /)
      atposMol(:,30) = (/   1.3599D0, -0.0138D0,-1.1667D0 /)
      atposMol(:,31) = (/  -3.4111D0,  2.4514D0,-0.0365D0 /)
      atposMol(:,32) = (/  -2.3959D0,  2.2441D0, 0.2666D0 /)
      atposMol(:,33) = (/  -3.7488D0,  3.3673D0, 0.4703D0 /)
      atposMol(:,34) = (/  -3.4494D0,  2.6027D0,-1.1248D0 /)
      atposMol(:,35) = (/  -4.8607D0, -1.3420D0, 0.0366D0 /)
      atposMol(:,36) = (/  -5.8759D0, -1.1347D0,-0.2666D0 /)
      atposMol(:,37) = (/  -4.8224D0, -1.4934D0, 1.1248D0 /)
      atposMol(:,38) = (/  -4.5230D0, -2.2579D0,-0.4703D0 /)
      !
      chargeMol = (/-0.180D0,  0.060D0,  0.060D0,  0.060D0, -0.120D0,  &
      &              0.060D0,  0.060D0, -0.120D0,  0.060D0,  0.060D0,  &
      &             -0.120D0,  0.060D0,  0.060D0, -0.060D0,  0.060D0,  &
      &             -0.060D0,  0.060D0, -0.120D0,  0.060D0,  0.060D0,  &
      &             -0.120D0,  0.060D0,  0.060D0, -0.120D0,  0.060D0,  &
      &              0.060D0, -0.180D0,  0.060D0,  0.060D0,  0.060D0,  &
      &             -0.180D0,  0.060D0,  0.060D0,  0.060D0, -0.180D0,  &
      &              0.060D0,  0.060D0,  0.060D0                       /)
      !                  Type,  atom 1, atom 2
      bondsMol(:, 1) = (/ 2, 2, 1 /)
      bondsMol(:, 2) = (/ 2, 3, 1 /)
      bondsMol(:, 3) = (/ 2, 4, 1 /)
      bondsMol(:, 4) = (/ 1, 5, 1 /)
      bondsMol(:, 5) = (/ 2, 6, 5 /)
      bondsMol(:, 6) = (/ 2, 7, 5 /)
      bondsMol(:, 7) = (/ 1, 8, 5 /)
      bondsMol(:, 8) = (/ 2, 9, 8 /)
      bondsMol(:, 9) = (/ 2,10, 8 /)
      bondsMol(:,10) = (/ 1,11, 8 /)
      bondsMol(:,11) = (/ 2,12,11 /)
      bondsMol(:,12) = (/ 2,13,11 /)
      bondsMol(:,13) = (/ 1,14,11 /)
      bondsMol(:,14) = (/ 2,15,14 /)
      bondsMol(:,15) = (/ 1,16,14 /)
      bondsMol(:,16) = (/ 2,17,16 /)
      bondsMol(:,17) = (/ 1,18,16 /)
      bondsMol(:,18) = (/ 2,19,18 /)
      bondsMol(:,19) = (/ 2,20,18 /)
      bondsMol(:,20) = (/ 1,21,18 /)
      bondsMol(:,21) = (/ 2,22,21 /)
      bondsMol(:,22) = (/ 2,23,21 /)
      bondsMol(:,23) = (/ 1,24,21 /)
      bondsMol(:,24) = (/ 2,25,24 /)
      bondsMol(:,25) = (/ 2,26,24 /)
      bondsMol(:,26) = (/ 1,27,24 /)
      bondsMol(:,27) = (/ 2,28,27 /)
      bondsMol(:,28) = (/ 2,29,27 /)
      bondsMol(:,29) = (/ 2,30,27 /)
      bondsMol(:,30) = (/ 1,31,14 /)
      bondsMol(:,31) = (/ 2,32,31 /)
      bondsMol(:,32) = (/ 2,33,31 /)
      bondsMol(:,33) = (/ 2,34,31 /)
      bondsMol(:,34) = (/ 1,35,16 /)
      bondsMol(:,35) = (/ 2,36,35 /)
      bondsMol(:,36) = (/ 2,37,35 /)
      bondsMol(:,37) = (/ 2,38,35 /)
      !
      !                   Type,  atom1,  atom2, atom3
      AnglsMol(:, 1) = (/ 3, 2, 1, 3 /)
      AnglsMol(:, 2) = (/ 3, 2, 1, 4 /)
      AnglsMol(:, 3) = (/ 2, 2, 1, 5 /)
      AnglsMol(:, 4) = (/ 3, 3, 1, 4 /)
      AnglsMol(:, 5) = (/ 2, 3, 1, 5 /)
      AnglsMol(:, 6) = (/ 2, 4, 1, 5 /)
      AnglsMol(:, 7) = (/ 2, 1, 5, 6 /)
      AnglsMol(:, 8) = (/ 2, 1, 5, 7 /)
      AnglsMol(:, 9) = (/ 1, 1, 5, 8 /)
      AnglsMol(:,10) = (/ 3, 6, 5, 7 /)
      AnglsMol(:,11) = (/ 2, 6, 5, 8 /)
      AnglsMol(:,12) = (/ 2, 7, 5, 8 /)
      AnglsMol(:,13) = (/ 2, 5, 8, 9 /)
      AnglsMol(:,14) = (/ 2, 5, 8,10 /)
      AnglsMol(:,15) = (/ 1, 5, 8,11 /)
      AnglsMol(:,16) = (/ 3, 9, 8,10 /)
      AnglsMol(:,17) = (/ 2, 9, 8,11 /)
      AnglsMol(:,18) = (/ 2,10, 8,11 /)
      AnglsMol(:,19) = (/ 2, 8,11,12 /)
      AnglsMol(:,20) = (/ 2, 8,11,13 /)
      AnglsMol(:,21) = (/ 1, 8,11,14 /)
      AnglsMol(:,22) = (/ 3,12,11,13 /)
      AnglsMol(:,23) = (/ 2,12,11,14 /)
      AnglsMol(:,24) = (/ 2,13,11,14 /)
      AnglsMol(:,25) = (/ 2,11,14,15 /)
      AnglsMol(:,26) = (/ 1,11,14,16 /)
      AnglsMol(:,27) = (/ 1,11,14,31 /)
      AnglsMol(:,28) = (/ 2,15,14,16 /)
      AnglsMol(:,29) = (/ 2,15,14,31 /)
      AnglsMol(:,30) = (/ 1,16,14,31 /)
      AnglsMol(:,31) = (/ 2,14,16,17 /)
      AnglsMol(:,32) = (/ 1,14,16,18 /)
      AnglsMol(:,33) = (/ 1,14,16,35 /)
      AnglsMol(:,34) = (/ 2,17,16,18 /)
      AnglsMol(:,35) = (/ 2,17,16,35 /)
      AnglsMol(:,36) = (/ 1,18,16,35 /)
      AnglsMol(:,37) = (/ 2,16,18,19 /)
      AnglsMol(:,38) = (/ 2,16,18,20 /)
      AnglsMol(:,39) = (/ 1,16,18,21 /)
      AnglsMol(:,40) = (/ 3,19,18,20 /)
      AnglsMol(:,41) = (/ 2,19,18,21 /)
      AnglsMol(:,42) = (/ 2,20,18,21 /)
      AnglsMol(:,43) = (/ 2,18,21,22 /)
      AnglsMol(:,44) = (/ 2,18,21,23 /)
      AnglsMol(:,45) = (/ 1,18,21,24 /)
      AnglsMol(:,46) = (/ 3,22,21,23 /)
      AnglsMol(:,47) = (/ 2,22,21,24 /)
      AnglsMol(:,48) = (/ 2,23,21,24 /)
      AnglsMol(:,49) = (/ 2,21,24,25 /)
      AnglsMol(:,50) = (/ 2,21,24,26 /)
      AnglsMol(:,51) = (/ 1,21,24,27 /)
      AnglsMol(:,52) = (/ 3,25,24,26 /)
      AnglsMol(:,53) = (/ 2,25,24,27 /)
      AnglsMol(:,54) = (/ 2,26,24,27 /)
      AnglsMol(:,55) = (/ 2,24,27,29 /)
      AnglsMol(:,56) = (/ 2,24,27,29 /)
      AnglsMol(:,57) = (/ 2,24,27,30 /)
      AnglsMol(:,58) = (/ 3,29,27,29 /)
      AnglsMol(:,59) = (/ 3,29,27,30 /)
      AnglsMol(:,60) = (/ 3,29,27,30 /)
      AnglsMol(:,61) = (/ 2,14,31,32 /)
      AnglsMol(:,62) = (/ 2,14,31,33 /)
      AnglsMol(:,63) = (/ 2,14,31,34 /)
      AnglsMol(:,64) = (/ 3,32,31,33 /)
      AnglsMol(:,65) = (/ 3,32,31,34 /)
      AnglsMol(:,66) = (/ 3,33,31,34 /)
      AnglsMol(:,67) = (/ 2,16,35,36 /)
      AnglsMol(:,68) = (/ 2,16,35,37 /)
      AnglsMol(:,69) = (/ 2,16,35,38 /)
      AnglsMol(:,70) = (/ 3,36,35,37 /)
      AnglsMol(:,71) = (/ 3,36,35,38 /)
      AnglsMol(:,72) = (/ 3,37,35,38 /)
      !                      Type, atom1, atom2, atom3, atom4
      DihrsMol(:, 1) = (/ 2,  6,  5,  1,  2 /)
      DihrsMol(:, 2) = (/ 2,  6,  5,  1,  3 /)
      DihrsMol(:, 3) = (/ 2,  6,  5,  1,  4 /)
      DihrsMol(:, 4) = (/ 2,  7,  5,  1,  2 /)
      DihrsMol(:, 5) = (/ 2,  7,  5,  1,  3 /)
      DihrsMol(:, 6) = (/ 2,  7,  5,  1,  4 /)
      DihrsMol(:, 7) = (/ 2,  8,  5,  1,  2 /)
      DihrsMol(:, 8) = (/ 2,  8,  5,  1,  3 /)
      DihrsMol(:, 9) = (/ 2,  8,  5,  1,  4 /)
      DihrsMol(:,10) = (/ 2,  9,  8,  5,  1 /)
      DihrsMol(:,11) = (/ 2,  9,  8,  5,  6 /)
      DihrsMol(:,12) = (/ 2,  9,  8,  5,  7 /)
      DihrsMol(:,13) = (/ 2, 10,  8,  5,  1 /)
      DihrsMol(:,14) = (/ 2, 10,  8,  5,  6 /)
      DihrsMol(:,15) = (/ 2, 10,  8,  5,  7 /)
      DihrsMol(:,16) = (/ 1, 11,  8,  5,  1 /)
      DihrsMol(:,17) = (/ 2, 11,  8,  5,  6 /)
      DihrsMol(:,18) = (/ 2, 11,  8,  5,  7 /)
      DihrsMol(:,19) = (/ 2, 12, 11,  8,  5 /)
      DihrsMol(:,20) = (/ 2, 12, 11,  8,  9 /)
      DihrsMol(:,21) = (/ 2, 12, 11,  8, 10 /)
      DihrsMol(:,22) = (/ 2, 13, 11,  8,  5 /)
      DihrsMol(:,23) = (/ 2, 13, 11,  8,  9 /)
      DihrsMol(:,24) = (/ 2, 13, 11,  8, 10 /)
      DihrsMol(:,25) = (/ 1, 14, 11,  8,  5 /)
      DihrsMol(:,26) = (/ 2, 14, 11,  8,  9 /)
      DihrsMol(:,27) = (/ 2, 14, 11,  8, 10 /)
      DihrsMol(:,28) = (/ 2, 15, 14, 11,  8 /)
      DihrsMol(:,29) = (/ 2, 15, 14, 11, 12 /)
      DihrsMol(:,30) = (/ 2, 15, 14, 11, 13 /)
      DihrsMol(:,31) = (/ 1, 16, 14, 11,  8 /)
      DihrsMol(:,32) = (/ 2, 16, 14, 11, 12 /)
      DihrsMol(:,33) = (/ 2, 16, 14, 11, 13 /)
      DihrsMol(:,34) = (/ 1, 31, 14, 11,  8 /)
      DihrsMol(:,35) = (/ 2, 31, 14, 11, 12 /)
      DihrsMol(:,36) = (/ 2, 31, 14, 11, 13 /)
      DihrsMol(:,37) = (/ 2, 17, 16, 14, 11 /)
      DihrsMol(:,38) = (/ 2, 17, 16, 14, 15 /)
      DihrsMol(:,39) = (/ 2, 17, 16, 14, 31 /)
      DihrsMol(:,40) = (/ 1, 18, 16, 14, 11 /)
      DihrsMol(:,41) = (/ 2, 18, 16, 14, 15 /)
      DihrsMol(:,42) = (/ 1, 18, 16, 14, 31 /)
      DihrsMol(:,43) = (/ 1, 35, 16, 14, 11 /)
      DihrsMol(:,44) = (/ 2, 35, 16, 14, 15 /)
      DihrsMol(:,45) = (/ 1, 35, 16, 14, 31 /)
      DihrsMol(:,46) = (/ 2, 19, 18, 16, 14 /)
      DihrsMol(:,47) = (/ 2, 19, 18, 16, 17 /)
      DihrsMol(:,48) = (/ 2, 19, 18, 16, 35 /)
      DihrsMol(:,49) = (/ 2, 20, 18, 16, 14 /)
      DihrsMol(:,50) = (/ 2, 20, 18, 16, 17 /)
      DihrsMol(:,51) = (/ 2, 20, 18, 16, 35 /)
      DihrsMol(:,52) = (/ 1, 21, 18, 16, 14 /)
      DihrsMol(:,53) = (/ 2, 21, 18, 16, 17 /)
      DihrsMol(:,54) = (/ 1, 21, 18, 16, 35 /)
      DihrsMol(:,55) = (/ 2, 22, 21, 18, 16 /)
      DihrsMol(:,56) = (/ 2, 22, 21, 18, 19 /)
      DihrsMol(:,57) = (/ 2, 22, 21, 18, 20 /)
      DihrsMol(:,58) = (/ 2, 23, 21, 18, 16 /)
      DihrsMol(:,59) = (/ 2, 23, 21, 18, 19 /)
      DihrsMol(:,60) = (/ 2, 23, 21, 18, 20 /)
      DihrsMol(:,61) = (/ 1, 24, 21, 18, 16 /)
      DihrsMol(:,62) = (/ 2, 24, 21, 18, 19 /)
      DihrsMol(:,63) = (/ 2, 24, 21, 18, 20 /)
      DihrsMol(:,64) = (/ 2, 25, 24, 21, 18 /)
      DihrsMol(:,65) = (/ 2, 25, 24, 21, 22 /)
      DihrsMol(:,66) = (/ 2, 25, 24, 21, 23 /)
      DihrsMol(:,67) = (/ 2, 26, 24, 21, 18 /)
      DihrsMol(:,68) = (/ 2, 26, 24, 21, 22 /)
      DihrsMol(:,69) = (/ 2, 26, 24, 21, 23 /)
      DihrsMol(:,70) = (/ 1, 27, 24, 21, 18 /)
      DihrsMol(:,71) = (/ 2, 27, 24, 21, 22 /)
      DihrsMol(:,72) = (/ 2, 27, 24, 21, 23 /)
      DihrsMol(:,73) = (/ 2, 29, 27, 24, 21 /)
      DihrsMol(:,74) = (/ 2, 29, 27, 24, 25 /)
      DihrsMol(:,75) = (/ 2, 29, 27, 24, 26 /)
      DihrsMol(:,76) = (/ 2, 29, 27, 24, 21 /)
      DihrsMol(:,77) = (/ 2, 29, 27, 24, 25 /)
      DihrsMol(:,78) = (/ 2, 29, 27, 24, 26 /)
      DihrsMol(:,79) = (/ 2, 30, 27, 24, 21 /)
      DihrsMol(:,80) = (/ 2, 30, 27, 24, 25 /)
      DihrsMol(:,81) = (/ 2, 30, 27, 24, 26 /)
      DihrsMol(:,82) = (/ 2, 32, 31, 14, 11 /)
      DihrsMol(:,83) = (/ 2, 32, 31, 14, 15 /)
      DihrsMol(:,84) = (/ 2, 32, 31, 14, 16 /)
      DihrsMol(:,85) = (/ 2, 33, 31, 14, 11 /)
      DihrsMol(:,86) = (/ 2, 33, 31, 14, 15 /)
      DihrsMol(:,87) = (/ 2, 33, 31, 14, 16 /)
      DihrsMol(:,88) = (/ 2, 34, 31, 14, 11 /)
      DihrsMol(:,89) = (/ 2, 34, 31, 14, 15 /)
      DihrsMol(:,90) = (/ 2, 34, 31, 14, 16 /)
      DihrsMol(:,91) = (/ 2, 36, 35, 16, 14 /)
      DihrsMol(:,92) = (/ 2, 36, 35, 16, 17 /)
      DihrsMol(:,93) = (/ 2, 36, 35, 16, 18 /)
      DihrsMol(:,94) = (/ 2, 37, 35, 16, 14 /)
      DihrsMol(:,95) = (/ 2, 37, 35, 16, 17 /)
      DihrsMol(:,96) = (/ 2, 37, 35, 16, 18 /)
      DihrsMol(:,97) = (/ 2, 38, 35, 16, 14 /)
      DihrsMol(:,98) = (/ 2, 38, 35, 16, 17 /)
      DihrsMol(:,99) = (/ 2, 38, 35, 16, 18 /)
      !
   case default
      info = 'Unsupport alkane, exit!'
      call error( subname, info, 1 )
   end select
   !
   molIDMol  = 1
   nMol      = 1
   molmass   = 0.D0
   do i = 1, natomMol
      molmass = molmass + EMass( attypMol(i) )
   enddo
   forall ( i=1:3 ) minbox(i) = maxval(atposMol(i,:)) - minval(atposMol(i,:))
   minlen = maxval( minbox )
   allocate( rotatedMol( 3, natomMol) )
   !
   subname = subnamep
return
end subroutine
